Ball mass - Ms = 2.1 t
Ball material - steel
Modulus of elasticity - Es = 206 GPa
Poisson′s ratio - νs = 0.3
Mass density - ρs = 7.85 tm3
Ball volume - Vs = Msρs = 2.1 t7.85 t/m3 = 0.268 m3 = 4/3πR³
Ball radius - Rs = 3 34·Vsπ = 400 mm
Height of bottom above the beam surface - H = 2 m
Structure type - simply supported beam
Beam length - L = 12 m
Material - reinforced concrete C20/25
Modulus of elasticity - E = 20 GPa
Poisson′s ratio - ν = 0.2
Shear modulus - G = E2· (1 + ν) = 20 GPa2· (1 + 0.2) = 8.33 GPa
Unit weight - γb = 25 kNm3
Cross section - rectangular with dimensions:
Width - b = 350 mm
Height - h = 650 mm
Area - A = b·h = 350 mm·650 mm = 2275 cm2
Second moment of area - I = b·h312 = 350 mm· (650 mm) 312 = 8009895833 mm4
Shear area - AQ = 56·A = 56·2275 cm2 = 1895.83 cm2
Self-weight - gb = A·γb = 2275 cm2·25 kN/m3 = 5.69 kN/m
Uniform load - q = 10 kNm
Gravity acceleration - g = 9.81 ms2
Uniform mass - m = gb + qg = 5.69 kN/m + 10 kN/m9.81 m/s2 = 1.6 t/m
The structure is reduced to a SDOF system for simplicity
Dynamically equivalent mass - Md = 2·Lπ·m = 2·12 m3.14·1.6 t/m = 12.22 t
Potential energy of the ball before dropping
Ep = Ms·g·H = 2.1 t·9.81 m/s2·2 m = 41.28 kJ
Kinetic energy immediately before the impact - Ek = Ms·v022
The velocity at the moment before the impact is determined by the energy conservation law Ek = Ep :
v0 = 2·EpMs = 2·41.28 kJ2.1 t = 6.26 m/s
Perfectly inelastic collision model is assumed.
Total mass after contact - Mtot = Ms + Md = 2.1 t + 12.22 t = 14.33 t
The velocity immediately after the contact is determined by the law of conservation of momentum:
v1 = v0·MsMtot = 6.26 m/s·2.1 t14.33 t = 0.92 m/s
Structural stiffness for a vertical force applied at the middle point of the span
K = 48·E·IL3 = 48·20 GPa·8009895833 mm4 (12 m) 3 = 4449.94 kN/m
Deflection due to uniform load
z0 = 5· (gb + q) ·L4384·E·I = 5· (5.69 kN/m + 10 kN/m) · (12 m) 4384·20 GPa·8009895833 mm4 = 26.44 mm
Static displacement - zst = Mtot·gK = 14.33 t·9.81 m/s24449.94 kN/m = 31.57 mm
Natural circular frequency - ω1 = KMtot = 4449.94 kN/m14.33 t = 17.62 s-1
Vibration period - T1 = 2·πω1 = 2·3.1417.62 s-1 = 0.356 s
Dynamic factor
μ = 1 + 1 + (v1·ω1g)2 = 1 + 1 + (0.92 m/s·17.62 s-19.81 m/s2)2 = 2.93
Dynamic displacement - zd = μ·zst = 2.93·31.57 mm = 92.58 mm
Dynamic force - Fd = μ·Ms·g = 2.93·2.1 t·9.81 m/s2 = 60.52 kN
(without self-weight and uniform load)
Simplified equation for the dynamic factor
μ1 = 1 + v1·ω1g = 1 + 0.92 m/s·17.62 s-19.81 m/s2 = 2.65
The difference will be smaller for greater heights.
Damped vibration is assumed with factor - ξ = 0.05
Vibration amplitude - A = zd − zst = 92.58 mm − 31.57 mm = 61.01 mm or
A = v1ω1 = 0.92 m/s17.62 s-1 = 52.21 mm
Theoretical equation of motion
y (t) = A·e-ξ·ω1·t·sin (ω1·t)
Solution by direct integration of the impulse load
Duration of impulse transmission for a beam with infinite mass [1]
τL = 2.94· 5 (1516·Ms·(1 − ν2E + 1 − νs2Es))2Rs·v0 = 2.94· 5 (1516·2.1 t·(1 − 0.2220 GPa + 1 − 0.32206 GPa))2400 mm·6.26 m/s = 3.93 ms
Duration of impulse transmission for a beam with finite mass [2]
τL = 2.94· 5 (1516·Ms·MdMs + Md·(1 − ν2E + 1 − νs2Es))2Rs·v0 = 2.94· 5 (1516·2.1 t·12.22 t2.1 t + 12.22 t·(1 − 0.2220 GPa + 1 − 0.32206 GPa))2400 mm·6.26 m/s = 3.69 ms
The above values correspond well to the experimental data in [3], where the recorded durations are of a similar magnitude.
The impulse force function will be determined by using the recommended expressions (9.20) - (9.22) in [1]
The coefficient of restitution for perfectly inelastic collision is - e. = 0
F (t) = Ms·v0· (1 + e.) ·π2·τL·sin(πτL·t)· (| t | ≤ τL)
Impulse load diagram
Maximal impulse load value - Fmax = F(τL2) = F(3.69 ms2) = 5613.66 kN
The equation of motion is expressed by the Duhamel′s integral
yD (t) = 1Mtot·ω1· min (t; τL) ∫0 ms F (τ) ·e-ξ·ω1· (t − τ) ·sin (ω1· (t − τ) ) dτ
Static displacement for the midpoint of the beam
y0 (t) = z0 + (zst − z0) ·{if t < T14: sin(2·π·tT1)
else: 1
Dynamic displacement for the midpoint of the beam
Number of intermediate joints - nJ = 11 (odd)
Length of one segment - Δx = LnJ + 1 = 12 m11 + 1 = 1 m
Coordinate of joint j - xJ (j) = Δx·j
Shear forces due to unit vertical load Fj = 1 at joint j
V1 (x; j) = {if x < xJ (j) : 1 − xJ (j) L
else: -xJ (j) L
Bending moments due to unit vertical load Fj = 1 at joint j
M1,max (j) = (xJ (j) L − 1)·xJ (j)
M1 (x; j) = M1,max (j) ·{if x < xJ (j) : xxJ (j)
else: L − xL − xJ (j)
Flexibility matrix
D (i; j) = ( L∫0 m M1 (x; i) ·M1 (x; j) dx)·1E·I + ( L∫0 m V1 (x; i) ·V1 (x; j) dx)·1G·AQ
D = 0.02160.03780.04890.05520.05740.0560.05140.04430.0350.0242⋯0.0124 0.03780.07040.0930.1060.1110.1090.10.08640.06850.0474⋯0.0242 0.04890.0930.1280.1490.1580.1550.1440.1240.09880.0685⋯0.035 0.05520.1060.1490.1790.1930.1930.180.1560.1240.0864⋯0.0443 0.05740.1110.1580.1930.2140.2170.2050.180.1440.1⋯0.0514 0.0560.1090.1550.1930.2170.2270.2170.1930.1550.109⋯0.056 0.05140.10.1440.180.2050.2170.2140.1930.1580.111⋯0.0574 0.04430.08640.1240.1560.180.1930.1930.1790.1490.106⋯0.0552 0.0350.06850.09880.1240.1440.1550.1580.1490.1280.093⋯0.0489 0.02420.04740.06850.08640.10.1090.1110.1060.0930.0704⋯0.0378 ⋮⋮⋮⋮⋮⋮⋮⋮⋮⋮⋱⋮0.01240.02420.0350.04430.05140.0560.05740.05520.04890.0378⋯0.0216 mm/kN
Mass matrix
dM.j = m·Δxt = 1.6 t/m·1 mt = 1.6
M = 1.6000000000⋯0 01.600000000⋯0 001.60000000⋯0 0001.6000000⋯0 00001.600000⋯0 000003.70000⋯0 0000001.6000⋯0 00000001.600⋯0 000000001.60⋯0 0000000001.6⋯0 ⋮⋮⋮⋮⋮⋮⋮⋮⋮⋮⋱⋮0000000000⋯1.6
Total mass of the structure - sum (⃗dM) = 19.7 t
Eigenvalues
Msq = √ M = 1.26000000000⋯0 01.2600000000⋯0 001.260000000⋯0 0001.26000000⋯0 00001.2600000⋯0 000001.920000⋯0 0000001.26000⋯0 00000001.2600⋯0 000000001.260⋯0 0000000001.26⋯0 ⋮⋮⋮⋮⋮⋮⋮⋮⋮⋮⋱⋮0000000000⋯1.26
C = copy (Msq·D·Msq; symmetric (nJ) ; 1; 1) = copy (Msq·D·Msq; symmetric (11) ; 1; 1) = 0.03450.06050.07810.08830.09180.1360.08220.07080.0560.0387⋯0.0198 0.06050.1130.1490.170.1780.2650.160.1380.110.0758⋯0.0387 0.07810.1490.2040.2380.2520.3780.230.1990.1580.11⋯0.056 0.08830.170.2380.2870.3090.4690.2870.250.1990.138⋯0.0708 0.09180.1780.2520.3090.3430.5290.3280.2870.230.16⋯0.0822 0.1360.2650.3780.4690.5290.8390.5290.4690.3780.265⋯0.136 0.08220.160.230.2870.3280.5290.3430.3090.2520.178⋯0.0918 0.07080.1380.1990.250.2870.4690.3090.2870.2380.17⋯0.0883 0.0560.110.1580.1990.230.3780.2520.2380.2040.149⋯0.0781 0.03870.07580.110.1380.160.2650.1780.170.1490.113⋯0.0605 ⋮⋮⋮⋮⋮⋮⋮⋮⋮⋮⋱⋮0.01980.03870.0560.07080.08220.1360.09180.08830.07810.0605⋯0.0345
⃗λ = reverse (last (eigenvals (C·10-3) ; 7) ) = [0.00261 0.000137 3.32×10-5 9.33×10-6 4.78×10-6 2.17×10-6 1.51×10-6]
Natural circular frequences - ⃗ω = 1⃗λ = [19.57 85.55 173.53 327.32 457.58 678.76 814.79] s⁻¹
Natural vibration frequences - ⃗f = ⃗ω2·π·Hz = ⃗ω2·3.14·Hz = [3.11 Hz 13.61 Hz 27.62 Hz 52.09 Hz 72.83 Hz 108.03 Hz 129.68 Hz]
Natural vibration periods - ⃗T = 1⃗f = [0.321 s 0.0734 s 0.0362 s 0.0192 s 0.0137 s 0.00926 s 0.00771 s]
Eigenvectors
Φ = inverse (Msq) ·extractcols (eigenvecs (C·10-3) ; range (nJ; nJ − 7 + 1; -1) ) = inverse (Msq) ·extractcols (eigenvecs (C·10-3) ; range (11; 11 − 7 + 1; -1) ) = 0.07520.1610.2170.280.3010.3230.309 0.1450.280.3190.280.186-9.74×10-13-0.113 0.2060.3230.2522.84×10-13-0.187-0.323-0.266 0.2530.280.0565-0.28-0.3091.15×10-130.216 0.2830.161-0.155-0.28-0.0270.3230.209 0.2930-0.251-7.46×10-150.225.01×10-14-0.193 0.283-0.161-0.1550.28-0.027-0.3230.209 0.253-0.280.05650.28-0.309-2.16×10-130.216 0.206-0.3230.252-2.77×10-13-0.1870.323-0.266 0.145-0.280.319-0.280.1861.09×10-12-0.113 ⋮⋮⋮⋮⋮⋮⋮0.0752-0.1610.217-0.280.301-0.3230.309
X = stack (matrix (1; 3) ; Φ; matrix (1; 3) ) = 0000000 0.07520.1610.2170.280.3010.3230.309 0.1450.280.3190.280.186-9.74×10-13-0.113 0.2060.3230.2522.84×10-13-0.187-0.323-0.266 0.2530.280.0565-0.28-0.3091.15×10-130.216 0.2830.161-0.155-0.28-0.0270.3230.209 0.2930-0.251-7.46×10-150.225.01×10-14-0.193 0.283-0.161-0.1550.28-0.027-0.3230.209 0.253-0.280.05650.28-0.309-2.16×10-130.216 0.206-0.3230.252-2.77×10-13-0.1870.323-0.266 ⋮⋮⋮⋮⋮⋮⋮0000000
Modal masses - ⃗mΦ = diag2vec (transp (Φ) ·M·Φ) ·t = [1 t 1 t 1 t 1 t 1 t 1 t 1 t]
Rayleigh damping model is assumed
β = 2·ξ⃗ω1 + ⃗ω2 = 2·0.0519.57 + 85.55 = 0.000951 , α = β·⃗ω1·⃗ω2 = 0.000951·19.57·85.55 = 1.59
ξ (ω) = α2·ω + β·ω2
Modal damping factors - ⃗ξΦ = ξ (⃗ω) = [0.05 0.05 0.0871 0.158 0.219 0.324 0.389]
Damped natural frequences
⃗ωD = ⃗ω· √ 1 − ⃗ξΦ^2·s-1 = [19.54 s-1 85.44 s-1 172.87 s-1 323.2 s-1 446.43 s-1 642.14 s-1 750.77 s-1]
Dynamic load vector
FΦ (i; t) = Φjm, i·F (t)
The equations of modal dynamic displacements are expressed by the Duhamel′s integral
yΦ (i; t) = 1⃗mΦ.i·ωD.i· min (t; τL) ∫0 ms FΦ (i; τ) ·e-⃗ξΦ.i·⃗ωi·s-1· (t − τ) ·sin (ωD.i· (t − τ) ) dτ
Joint displacements
yJ (j; t) = 7∑i= 1Φj, i·yΦ (i; t)
Comparison of time history records at midpoint for SDOF and MDOF systems
Time history records for separate joints
Beam deflections at the first five time steps of Δt = 17.7 ms
Animation of beam response (slowed down)
[1] Harris C. M., Piersol A.G., HARRIS’ SHOCK AND VIBRATION HANDBOOK, Fifth Edition, McGraw-Hill 2002, ISBN 0-07-137081-1
[2] Qing Peng, Xiaoming Liu, Yueguang Wei, Elastic impact of sphere on large plate, Journal of the Mechanics and Physics of Solids, Volume 156, 2021, 104604, ISSN 0022 - 5096, https://doi.org/10.1016/j.jmps.2021.104604
[3] Hong Hao and Thong M. Pham, Performance of RC Beams with or without FRP Strengthening Subjected to Impact Loading, Proceedings of the 2nd World Congress on Civil, Structural, and Environmental Engineering (CSEE’17) Barcelona, Spain – April 3– 4, 2017 ISSN:2371 - 5294 DOI:10.11159/icsenm17.1
[4] Gugan, D. “Inelastic collision and the Hertz theory of impact.” American Journal of Physics 68 (2000): 920-924., http://www.oxfordcroquet.com/tech/gugan/index.asp
